{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# State space representation\n",
    "\n",
    "The \"standard\" or most commonly used state space representation is \n",
    "\n",
    "\\begin{align} \\dot{x} &= Ax + Bu \\\\ y &= Cx + Du \\end{align}\n",
    "\n",
    "Take note that Seborg uses a slightly different version:\n",
    "\n",
    "\\begin{align} \\dot{x} &= Ax + Bu + Ed\\\\ y &= Cx \\end{align}\n",
    "\n",
    "This second version can not represent pure gain systems as it effectively assumes $D=0$. It is also possible to stack $u$ and $d$ from the bottom form into one input vector, so the $E$ matrix really doesn't add much. As you may infer, I prefer the top version and it is also the version used by most libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy\n",
    "import numpy.linalg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Converting between state space and transfer function forms\n",
    "\n",
    "There is good support in various libraries for converting systems with numeric coefficients between transfer function and state space representation."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Scipy.signal\n",
    "\n",
    "The `scipy.signal` library handles conversion between transfer function coefficients and state space matrices easily. Note that `scipi.signal` only handles SISO transfer functions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import scipy.signal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "G = scipy.signal.lti(1, [1, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "TransferFunctionContinuous(\n",
       "array([1.]),\n",
       "array([1., 1.]),\n",
       "dt: None\n",
       ")"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This object allows us to access the numerator and denominator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([1.]), array([1., 1.]))"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G.num, G.den"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To convert to state space, we can use the `.to_ss()` method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "Gss = G.to_ss()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[-1.]]), array([[1.]]), array([[1.]]), array([[0.]]))"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gss.A, Gss.B, Gss.C, Gss.D"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can build another object using the state space matrices instead of the Laplace form"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "StateSpaceContinuous(\n",
       "array([[-1.]]),\n",
       "array([[1.]]),\n",
       "array([[1.]]),\n",
       "array([[0.]]),\n",
       "dt: None\n",
       ")"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G2ss = scipy.signal.lti(Gss.A, Gss.B, Gss.C, Gss.D)\n",
    "G2ss"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can convert to transfer function form using `.to_tf()` (there is a small warning about bad coefficients, but the answer is reliable)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/homebrew/anaconda3/envs/dynamicscontrol/lib/python3.11/site-packages/scipy/signal/_filter_design.py:1746: BadCoefficients: Badly conditioned filter coefficients (numerator): the results may be meaningless\n",
      "  warnings.warn(\"Badly conditioned filter coefficients (numerator): the \"\n"
     ]
    }
   ],
   "source": [
    "G2 = G2ss.to_tf()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now access the numerator and denominator again:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([1.]), array([1., 1.]))"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G2.num, G2.den"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Instead of building objects we can also use the functions in `scipy.signal.lti_conversion`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/var/folders/0n/nkm9qrn50h5gj7_5qc_tqzvc0000gn/T/ipykernel_60648/1776210230.py:1: DeprecationWarning: Please use `tf2ss` from the `scipy.signal` namespace, the `scipy.signal.lti_conversion` namespace is deprecated.\n",
      "  scipy.signal.lti_conversion.tf2ss(1, [1, 1])\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(array([[-1.]]), array([[1.]]), array([[1.]]), array([[0.]]))"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "scipy.signal.lti_conversion.tf2ss(1, [1, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/var/folders/0n/nkm9qrn50h5gj7_5qc_tqzvc0000gn/T/ipykernel_60648/778838862.py:1: DeprecationWarning: Please use `ss2tf` from the `scipy.signal` namespace, the `scipy.signal.lti_conversion` namespace is deprecated.\n",
      "  scipy.signal.lti_conversion.ss2tf(-1, 1, 1, 0)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(array([[0., 1.]]), array([1., 1.]))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "scipy.signal.lti_conversion.ss2tf(-1, 1, 1, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Control library\n",
    "\n",
    "The control library (at least from version 0.8.0) does a good job with these conversions as well."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "import control"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\frac{1}{s + 1}$$"
      ],
      "text/plain": [
       "TransferFunction(array([1]), array([1, 1]))"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gtf = control.tf([1], [1, 1])\n",
    "Gtf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the control library we convert the system using `ss` (short for state space) to get a State Space representation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\n",
       "\\left(\\begin{array}{rll|rll}\n",
       "-1\\phantom{.}&\\hspace{-1em}&\\hspace{-1em}\\phantom{\\cdot}&1\\phantom{.}&\\hspace{-1em}&\\hspace{-1em}\\phantom{\\cdot}\\\\\n",
       "\\hline\n",
       "1\\phantom{.}&\\hspace{-1em}&\\hspace{-1em}\\phantom{\\cdot}&0\\phantom{.}&\\hspace{-1em}&\\hspace{-1em}\\phantom{\\cdot}\\\\\n",
       "\\end{array}\\right)\n",
       "$$"
      ],
      "text/plain": [
       "<LinearIOSystem:sys[2]:['u[0]']->['y[0]']>"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gss = control.ss(Gtf)\n",
    "Gss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-1.]])"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gss.A"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Symbolic conversion"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It is easy to convert state space models to transfer functions since the Laplace transform is a linear operator:\n",
    "\n",
    "$$ \\dot{x} = Ax + Bu \\quad \\therefore \\quad sX(s) = AX(s) + BU(s) \\quad X(s) = (sI - A)^{-1}BU(s)$$\n",
    "$$ y = Cx + Du \\quad \\therefore \\quad Y(s) = CX(s) + DU(s) \\quad Y(s) = \\underbrace{(C(sI - A)^{-1}B + D)}_{G(s)}U(s)$$\n",
    "\n",
    "This conversion is handled for symbolic matrices by `tbcontrol.symbolic.ss2tf`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sympy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "import tbcontrol\n",
    "tbcontrol.expectversion('0.1.8')\n",
    "import tbcontrol.symbolic"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = sympy.symbols('s')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "A, B, C, D = [sympy.Matrix(m) for m in [G2ss.A, G2ss.B, G2ss.C, G2ss.D]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(Matrix([[-1.0]]), Matrix([[1.0]]), Matrix([[1.0]]), Matrix([[0]]))"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A, B, C, D"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}\\frac{1.0}{s + 1.0}\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([[1.0/(s + 1.0)]])"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G = tbcontrol.symbolic.ss2tf(A, B, C, D, s)\n",
    "G"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that `ss2tf` returns a sympy Matrix. To get the SISO result, we need to index into the matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{1.0}{s + 1.0}$"
      ],
      "text/plain": [
       "1.0/(s + 1.0)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G[0, 0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice that the roots of the characteristic function correspond with the eigenvalues of the A matrix. The numerator and denominator of control transfer functions are stored as lists of lists to accomodate MIMO systems."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([-1.+0.j])"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gtf.poles()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([-1.])"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "numpy.roots(Gtf.den[0][0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "EigResult(eigenvalues=array([-1.]), eigenvectors=array([[1.]]))"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "numpy.linalg.eig(Gss.A)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}